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We present a new covariant, gauge-invariant formalism describing linear metric perturbation fields 
on any spherically symmetric background in general relativity. The advantage of this formalism re- 
lies in the fact that it does not require a decomposition of the perturbations into spherical tensor 
harmonics. Furthermore, it does not assume the background to be vacuum, nor does it require its 
staticity. In the particular case of vacuum perturbations, we derive two master equations describ- 
ing the propagation of arbitrary linear gravitational waves on a Schwarzschild black hole. When 
decomposed into spherical harmonics, they reduce to covariant generalizations of the well-known 
Regge- Wheeler and Zerilli equations. Next, we discuss the general case where the metric per- 
turbations are coupled to matter fields and derive a new constrained wave system describing the 
' propagation of three gauge-invariant scalars from which the complete metric perturbations can be 

, reconstructed. We apply our formalism to the Einstein-Euler system, dividing the fluid perturba- 

tions into two parts. The first part, which decouples from the metric perturbations, obeys simple 
advection equations along the background flow and describes the propagation of the entropy and 
the vorticity. The second part describes a perturbed potential flow, and together with the metric 
perturbations it forms a closed wave system. 
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I. INTRODUCTION 



O 

The purpose of this work is to develop a gauge-invariant perturbation formalism describing the propagation of 
5^ | linearized gravitational and matter fields on an arbitrary spherically symmetric background configuration. Examples 
of such configurations include nonrotating black holes and stars, spherical matter distributions undergoing complete 
gravitational collapse, or expanding shells of matter. Therefore, it is clear that a perturbation formalism for such 
backgrounds possesses a wide range of applications, covering the stability analysis of compact objects, the description 
of quasi-normal oscillations of these objects, the computation of linearized gravitational waves produced by a small star 
moving in the field of a nonrotating black hole or by another black hole in the close limit approximation, the stability 
C*") ■ analysis of Cauchy horizons in gravitational collapse models, and the linear evolution of structures in cosmological 
f"*** | models. 1 

Linear gravitational fluctuations on a Schwarzschild black hole were first studied by Regge and Wheeler 0] and 
0^ ' Zerilli Q, who decomposed the perturbations into spherical tensor harmonics, analyzed their behavior with respect 
to infinitesimal coordinate transformations, and by fixing the gauge, obtained a family of master equations describing 
odd- and even-parity perturbations with a given angular momentum number i > 2. These master equations, known 
as the Regge- Wheeler and the Zerilli equations, respectively, describe the propagation of the true dynamical degrees 
of freedom of the linear theory, completely eliminating gauge and constraint modes. Therefore, they are ideally 
suited for understanding the behavior of linear perturbations of the Schwarzschild black hole. From the fact that 
these equations can be written as a 1 + 1-dimensional wave equation with a positive, time-independent potential, the 
existence of exponentially growing modes is immediately ruled out and stability in this sense follows. With a little 
bit more of work, and based on the techniques described in Refs. Q and Q, one can also prove that the solutions 
to the Regge- Wheeler and Zerilli equations belonging to sufficiently regular initial data on a spacelike slice remain 
uniformly bounded outside the black hole. 

The derivation of the Regge- Wheeler and Zerilli master equations was later clarified by the work of Moncricf 6] 
who analyzed the gravitational perturbations from a Hamiltonian point of view and casted the Regge- Wheeler and 
Zerilli equations into gauge-invariant form. Moncrief also obtained master equations describing linear perturbations of 
spherically symmetric fluid stars Q and linear gravito-electromagnetic perturbations of a Reissner-Nordstrom black 
hole @t13i an d showed that the latter is linearly stable in Einstein-Maxwell theory. 2 

Later, Gerlach and Sengupta [l2j provided a covariant description of the gauge-invariant perturbation approach, 
where instead of assuming a Schwarzschild background written in the usual Schwarzschild coordinates (t, r) and 
foliated by the static t — const slices, they only assumed the background to be spherically symmetric. Based on 



1 See Ref. Q for a review on perturbation theory for Schwarzschild black holes and references to some of the applications we mention. 

2 Interestingly, the magnetically charged Reissner-Nordstrom black hole is linearly unstable in Einstein- Yang-Mills-Higgs theory [ill . 
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the natural 2 + 2 form of the background induced by the metric two-spheres and the two-dimensional orbit manifold 
M = M/SO(3) orthogonal to them, and based on the decomposition into spherical tensor harmonics, Gerlach and 
Sengupta introduced a complete set of gauge-invariant tensor fields on the manifold M. This covariant approach has 
several advantages. Besides the possibility of describing perturbations of dynamical spacctimes, such as collapsing 
spherical stars, it also allows to describe the propagation of linearized gravitational waves on a Schwarzschild black 
hole in local coordinates which arc regular at the event horizon, such as Kruskal- or ingoing Eddington-Finkclstein 
coordinates. Based on their approach, Gerlach and Sengupta obtained the covariant form of the Regge- Wheeler master 
equation, valid in any coordinate system compatible with the 2 + 2 form of the background. The covariant form of the 
Zerilli master equation was derived in [l3| (see also [l4[ ) and a generalization including source terms with applications 
to the problem of calculating gravitational waves produced by the motion of a small star moving around a black hole 
was given in Rcf. [l5|. The covariant, gauge-invariant approach has also been applied to the derivation of master 
equations describing linear perturbations of black holes in general relativity coupled to a nonlinear electromagnetic 
theory, see Ref . (l6j where sufficient conditions for the linear stability of such holes are also given. Further developments 
of the covariant gauge-invariant formalism include the coupling to a perfect fluid (l7j and the generalization to second 
and higher-order perturbation theory p^ - [2l| . 

In a somewhat different development, linear perturbations of a Schwarzschild black hole were also analyzed in [22[ 
based on the 3+1 formulation. The main advantage of this work is that, unlike previous approaches, the decomposition 
into spherical tensor harmonics is not needed, simplifying the derivation of the equations. Instead, two canonical pairs 
of gauge- invariant scalars are constructed, describing axial and polar perturbations in the " mono-dipolc-frec" sector 
(i.e. those fields corresponding to £ > 2 in the decompositions into spherical harmonics). In the polar case, these 
gauge- invariant scalars involve quasilocal operators which are local on M, but non-local on the two-spheres S 2 . Based 
on these quantities, two scalar wave equations are derived which reduce to the covariant forms of the Regge- Wheeler 
and Zerilli master equations when decomposed into spherical harmonics. While the two scalar equations are presented 
in the covariant description, their derivation is based on the Schwarzschild coordinate patch. 

Yet a different approach to gravitational perturbation theory which does not require the background to be spherically 
symmetric but assumes instead that it is static casts the perturbation equations into a wave equation for the linearized 
extrinsic curvature tensor [ij], [23Tj25j . This curvature-based approach has turned out to be useful for establishing the 
linear stability of certain Einstein- Yang-Mills black holes with a negative cosmological constant [111 [27} . 

In this work, we combine the covariant, gauge-invariant approach of fl2j with the quasilocal method in [22j and 
present a covariant, gauge-invariant perturbation formalism for an arbitrary spherically symmetric background without 
performing the expansion into spherical tensor harmonics. To this purpose, we first review the relevant background 
equations describing the most general spherically symmetric spacctimc in Sec. [TTJ Then, as a warmup, we derive in 
Sec. lIIII master equations describing the propagation of scalar and electromagnetic test fields on such spacetimes. These 
master equations have the form of a wave equation on M with an effective potential and act on an angular-dependent 
scalar field on M which is gauge-invariant in the electromagnetic case. Next, in Sec. IIVI we discuss linear metric 
perturbation and first review their behavior under infinitesimal coordinate transformations. Based on the quasilocal 
approach and a decomposition of tensor fields on the sphere in terms of scalars which is discussed in the appendix, 
we construct a set of gauge-invariant, angular-dependent tensor fields on M which behave as scalars under rotations 
of the two-spheres. When performing an expansion into spherical harmonics, they reduce to the gauge-invariant 
tensors introduced by Gerlach and Sengupta. Then, we derive the expressions for the linearized Riemann curvature, 
Ricci and Einstein tensors, and obtain the perturbed Einstein equations. Next, in Sec. [V]we consider the vacuum 
case and derive the gauge-invariant master wave equations which reduce to the covariant Regge- Wheeler and Zerilli 
equations when an expansion in spherical harmonics is performed. We also show that the Regge- Wheeler equation 
can naturally be obtained in both the odd- and the even-parity sectors and discuss how the metric perturbations can 
be reconstructed from the scalar potentials satisfying the master equations. 

Next, in Sec. I VII we consider the coupling of the metric fields to arbitrary matter fields. In the odd-parity sector, 
we discuss the generalization of the Regge- Wheeler equation and state assumptions under which it yields a master 
equation. In the even-parity sector, we do not derive master equations for the metric fields; instead, we derive a 
constrained wave system for two gauge-invariant scalars from which the metric perturbations can be reconstructed 
under certain assumptions on the matter fields which should be reasonable. In the vacuum case, the two wave equations 
decouple from each other and are related to the radial parts of the Teukolsky equations [28| for the Weyl scalars 
and $4. In fact, as we explain, our equations are equivalent to the ones obtained by Bardeen and Press [29(, where 
is called ^2 and M/4 is called ^-2- Finally, in Sec. IVIIl we apply our formalism to the perturbations of self-gravitating 
spherical fluid configurations. Focussing on the fluid perturbations first, we decompose them into two parts, where the 
first part determines the perturbed vorticity and entropy and the second part describes a perturbed potential flow. 
The propagation of the first part is described by simple advection equations along the background flow. It decouples 
completely from the second part and the metric perturbations. Therefore, it can be solved separately. The second 
part couples to the metric perturbations, and showing that the assumptions of the previous section are satisfied, we 
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derive a system of effective wave equations on the orbit manifold M, describing the propagation of the remaining 
part of the even-parity perturbations. This system should be useful for analytic and numerical investigations of the 
aforementioned problems. A summary and conclusions are given in Sec. I Villi 

We use the signature convention (— , +, +, +) for the metric and choose units for which c = 1. 



II. BACKGROUND EQUATIONS 

A spherically symmetric spacetime (M, g) can be written as the product of a two-dimensional pscudo-Ricmannian 
manifold (M , g) with the two-sphere (5' 2 ,g), 

M = Mx S 2 , g = g ab dx a dx b + r 2 g AB dx A dx B , (1) 

where r is a strictly positive function on M, and x a ,x b , . . . and x A ,x B , . . . denote local coordinates on M and S 2 , 
respectively. The geometric interpretation of the function r is the following: let p £ M and consider the two-sphere 
S p := {p} x S 2 C M. Let A{p) = \S P \ denote the area of S p , computed from the induced metric on S p . Therefore, 
A(p) = A-Kr(p) 2 and the function r is defined geometrically as 



and is called the areal radius. For the following, we assume M to be orientable such that it has the volume form 
£ab '■= \ZW\ e ab with \g\ :— I det(g a b) | the absolute value of the determinant of g a b and eoo = en = 0, £oi = — eio = 1- 
We also introduce the differential r a := V Q r and the covariant Hessian r a b '■= V a Vbr of r with respect to the covariant 
derivative V associated to the two-metric g. 

The Christoffel symbols corresponding to the metric in Eq. (jTJ) arc 

T d ab = f d ab , T d aB = 0, T d AB = -rr d g AB , (2a) 

^ D ab = 0, r fl „B = — S D B , T D AB~^AB, (2b) 

r 

where here and in the following, quantities with a tilde and a hat refer to the manifolds (M , g) and (S, g), respectively. 
From this, one finds the following expressions for the curvature tensor, 

-Jbd - 9ad9bc), (3a) 

(3b) 

rr ac g BDl (3c) 

(3d) 
(3c) 

j9cb), (3f) 

where k denotes the Gauss curvature of (M, g) and where N := g(dr,dr) = r a r a . Contracting, one obtains the Ricci 
tensor, 

!— , (4a) 
r 

(4b) 

AB j (4c) 



Rabcd 


= R, 


RaBcd 


= o, 


RaBcD 


= — \ 


RabCD 


= o, 


RaBCD 


= o, 


Rabcd 


= r 2 



Rab 


= kgab 


RaB 


= o, 


Rab 


= ('- 



where Ar := g r a b = V a V a ?'. From this, one finally obtains the Einstein tensor, 



Gab = -- (r ab ) tf ~ A,g ab (1-N-rAr), (5a) 



2, , ft 1_ 

2 ; 

G aB = 0, (5b) 
..2' 



Gab = rAr - kr 2 g AB , (5c) 
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where (r a b)*^ '■= T ab — \g a bg cd Tcd denotes the trace-free part of r ab . This particular structure of the Einstein tensor im- 
plies that the stress-energy tensor T must satisfy the conditions T aB = and T AB proportional to g AB . Summarizing, 
Einstein's field equations in spherical symmetry consist of 

~(r ab ) tf = kTH (6a) 
~(l-N- rAr) = K g ab T ab , (6b) 
2 (rAr - I-r 2 ) = K g AB T AB , (6c) 

where k := 8itGn with Newton's constant Gn- 
The Bianchi identities imply that 

= V M G^ b = lv a (r 2 G a h ) - ^g AB G AB . 

If the stress-energy tensor is divergence-free, the same equation holds for T, and hence, 

jg AB (G AB - kT ab ) = V a [r 2 (G a b - KT a b )} . 

This shows that Eq. (|6c| follows from Eqs. (|6a|6b|) provided r b ^ and the divergence-free equation for T holds. 
Therefore, it is sufficient to solve Eqs. (|6al6bj) . 



III. SCALAR AND ELECTROMAGNETIC FIELDS PROPAGATING ON A SPHERICALLY 

SYMMETRIC BACKGROUND 

In this section we derive the master equations describing scalar and electromagnetic fields propagating on an 
arbitrary spherically symmetric spacetime. 



A. Scalar field propagation 

The propagation of a scalar field on a fixed spacetime background (M , g) is described by Klein-Gordon equation 

(□ g + M 2 )0 = O, (7) 

where is a real scalar field on M with mass /i (in units for which h = c = 1), and D g := — V M V M is the covariant 
d'Alcmbcrtian on (M, g). Assuming that (M, g) is spherically symmetric and applying the general formula 



1 d_ ( /TT uu d( t> 



= ( vW"^ ) . \g\ ■■= I det( 5 ^) 



to a local coordinate patch, we can rewrite Eq. ([7} in its 2 + 2 form, 

_ J_v Vv a 0) + (~ + ^ = 0, 

where A = g AB V b is the Laplacian on the round sphere (5 2 ,g). Finally, we introduce the rescaled scalar field 
ip := rcj) and use the identity 

-V Q (r 2 V a 0) = rDiP + (Ar)ip, 
with □ := —A = — V a V a the covariant d'Alcmbcrtian on (M,g). This yields 

/ A Ar \ 

Di/)+ — j + — +m 2 V' = 0. (8) 
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This has the form of a wave equation on (M, g) with effective potential 

A Ar 2 
V : =-^2+ — +f 1 ■ 

Here, the first term is an operator on the sphere S 2 . If ijj is decomposed into spherical harmonics, it becomes the 
multiplicative operator £(£+!) /r 2 with I the angular momentum number, which represents the usual centrifugal term. 
The second term is a curvature correction term which can be further simplified using the Einstein equation (|6b[) . For 
a Schwarzschild spacetime of mass m, for instance, for which N = 1 — 2m/r, we have r _1 Ar = 2m/r 3 . The third 
term in the potential, [i 2 , is just inherited from the corresponding term of the original Klein-Gordon equation (J7|). 

As we will see below, the effective equations describing electromagnetic and linearized fluctuations on spherically 
symmetric spacetimes have a very similar form than Eq. ([8]) . 



B. Electromagnetic propagation 

Next, we consider Maxwell's equations 

- VF^ = J,, iV = V„A, - (9) 

on a spherically symmetric background (M, g). Here, A = A^dx^ is the electromagnetic potential one-form, F := 
dA = 2 _1 F^dx^ A dx v the corresponding Faraday tensor, and J = J^dx^ is the electric four-current density. 
According to the 2 + 2 form of the background, we may decompose the electromagnetic potential as 

A = a a dx a + j3 B dx B , 

where the quantities a a and /3b both depend on the coordinates (x a ,x B ). Notice that with respect to a rotation 
on S 2 , Pa transforms like the components of a one- form while a a transforms as a scalar. With respect to a gauge- 
transformation A i— > A + dt; parametrized by a function £ on M, we have 

a a h> a a + Vat, Pb^Pb + Vb£- (10) 

Quantities which are invariant with respect to these transformations can be constructed based on the following 
decomposition for a one- form uj = u>Adx A in terms of two scalar fields / and g on the two-sphere S 2 : 

cuA = V A f + iA B V B g 7 (11) 

where sab denotes the volume form on S 2 . See the appendix for a proof and further discussion on this decomposition. 
Here, we note that the two terms on the right-hand side of Eq. (|TT|) are mutually orthogonal with respect to the natural 
L 2 scalar product, 

:= J ' g AB u A r}ByAF\d 2x 

s 2 

for two one- forms u) and r/ on S 2 . Therefore, the decomposition (1111) is unique. However, the functions / and g 
themselves are only unique up to an additive constant. In the following, we fix this constant by requiring / and g to 
have zero mean values over (S 2 ,g). 

The decomposition allows us to represent 

Pb = Vb/U + e B C Vcv 

in terms of two scalar fields /i and v. Since /3b depends not only on the angular variables, but also on the radial ones, 
the fields fi and v depend on both x a and x B . However, unlike (3b , they transform like scalar fields under rotations of 
the two-spheres. Similarly, the one-form a a dx a depends on the angular coordinates x B , but it transforms as a scalar 
field under rotations. As mentioned above, from now on we assume that fi and v have zero mean values over (S" 2 ,^). 
This implies that their decompositions into spherical harmonics have a vanishing monopole term. 

Working in the monopole- free space, the gauge transformations (|10p imply the following transformations for a a , fi 
and v. 



a a >-» oi a + V a £, /x h-> ii + £, v H- v. 



(12) 
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Therefore, the following quantities are gauge-invariant, i.e. invariant with respect to the transformations (|10|) : 

ai mv *> := a a - V a /i, V (mv) := i/. (13) 

The Faraday tensor F is gauge- invariant, and hence it must be possible to express it solely in terms of aa and 
p(mv) _ short calculation reveals that 

^ab — V a a fc — Vf,a^ i* a s — — t Ba — — Vfla> + £b VcV„i/ 1 b ab — —£AB<±v y 

Applying the general formula = —F ufl ) 

to a coordinate patch, we obtain the following Maxwell equations in their 2 + 2 form, 

= r 2 J b , (14a) 
= 3, (14b) 
= k, (14c) 

where we have decomposed the angular components of the four-current according to Jb = VbJ + £s c 'Vcfc. Eq. (|14c[) 
decouples from the remaining two equations, and has the form of an inhomogeneous wave equation on (M, g) with 
effective potential V = — r~ 2 A. This master equation describes the propagation of the odd-parity (or axial) part of 
the electromagnetic field. 

In order to derive a master equation for the even-parity (or polar) part of the field, we first introduce the scalar 

in terms of which F a t = —r~ 2 i a b& and Eq. (|14a[) can be rewritten as 

e ab V Q $ - Aa { b mv) = r 2 J b . (15) 
Applying the two-dimensional curl operator e ab W a on both sides of this equation we obtain 

□$--^A$ = -^e ab V a (r 2 J fc ), (16) 

which is an inhomogeneous wave equation on (M, g) of exactly the same form as the odd-parity equation (|14c[) . It 
is a master equation describing the propagation of the even-parity part of the field: solving Eq. (|16[) gives $, from 

which the gauge-invariant quantity (x£ nv ^ can be reconstructed from Eq. (|15p since the Laplacian A is invertible on 
the monopole-free space. This quantity automatically satisfies Eq. (|14b|) since by Eq. (|15[) and the continuity equation 
V^J^ = we have — A(V b a i b " l "' > ) = V b (? ,2 Jb) = — Aj. Therefore, we can reconstruct the even-parity part of the 
Faraday tensor, and all the Maxwell equations are satisfied. 

IV. GRAVITATIONAL PERTURBATION EQUATIONS 

In this section we describe linearized metric perturbations on a spherically symmetric background. We start by 
analyzing their behavior under an infinitesimal coordinate change, and based on ideas introduced in the previous 
section, construct a full set of gauge-invariant, angular-dependent tensor fields on M. Similar invariants are also 
constructed from the perturbed stress-energy tensor. Then, we derive the linearized Einstein equations and write 
them in terms of these gauge-invariant quantities. 

We consider a smooth perturbation of the spherically symmetric metric g, that is, a smooth, one-parameter family 3 
of metrics, g(A), on M such that g(0) = g. To first order, the deviation from the background metric g is described 



3 It is also possible to perturb M , considering a one-parameter family of spacetimes (M(A), g(A)). In this case, the manifolds M(A) need 
to be identified with each other by an appropriate diffcomorphism, which leads to an alternative, but equivalent point of view. For a 
discussion which is based on this approach and analyzes the behavior of linear and nonlinear perturbations under diffcomorphisms, see 
Ref. [3dl for instance. 



V7 (inv) 



Ai 



(inv) 



ya a (in«) 
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by the variation 



Sg := 



d 

dX* 



(A) 



A=0 



whose geometric interpretation is the tangent vector to the curve g(A) at the point g = g(0) in configuration space. 
With respect to a one-parameter family of diffcomorphisms, <f>(X) : M — > M, the metrics g(A) transform according to 



g(A)^0(A)*g(A), 



(17) 



the star denoting pull-back. We may split 4>{\) = xp(X)<p, where <p := (f)(0) : M -+ M is the diffcomorphism to zeroth 
order in A and where : M -+ M satisfies ^(0) = id. Differentiating Eq. ()17|) with respect to A, we then obtain 



<5g H' if*(5g + £ x g), 
where the vector field X is defined as the variation of tp, 



(18) 



* = ^(A)(p) 



A=0 



For the following, we will restrict ourselves to " background" diffeomorpisms <p which leave the structure M = M x S 2 
of the manifold invariant. The next step in the discussion of our perturbation formalism is to construct quantities 
from <5g which arc invariant with respect to the infinitesimal coordinate transformations <5g t— > Sg + ,£ x g. 



A. 2 + 2 split and gauge invar iance 



For the following, it is convenient to split the metric perturbations 5g in accordance with the 2 + 2 form of the 
background metric: 

H a b ■= Sg ab , Q aB ■= Sg a B , K A b ■= r~ 2 5gAB , 

where the quantities H a b, Q a B and Kab depend on the coordinates (x a , x B ). Similarly, we may split the vector field 
X in Eq. (fTH)) according to 

i a :=X a , if:=X B . 
With respect to the infinitesimal coordinate transformation generated by X we then have 

H ab h+ H ab + 2V (a e fc) , (19a) 

QaB ^ QaB + Vs^a + r 2 V a f te , (19b) 

2 

Kab i-> Kab + - <}ab r a £ a + 2\7 {a Vb) , (19c) 

with f] A ■= 9abV B 

In order to construct quantities which are invariant with respect to these transformations, we generalize the ideas 
presented in the previous section and decompose any tensor field on 5* 2 in terms of scalars. For one-forms oj on S* 2 
we recall the decomposition in Eq. (|11|) above. For a symmetric, trace-free tensor field r = TAsdx A dx B on S 2 we use 
the following decomposition proved in the appendix, 

tab = (VaVbF) tf + e (A C V B) VcG, (20) 

where F and G are scalar fields on S 2 and the super-index TF denotes the trace-free part with respect to g. This 
decomposition is also orthogonal with respect to the natural L 2 scalar product on S 2 and therefore, it is unique. 
However, the functions F and G are only unique up to the addition of a monopole or dipole term (that is, a function 
L on S 2 satisfying (V aV ' bL) tf = 0). We fix this freedom by working on the mono-dipole-free space defined as the 
orthogonal complement of the monopole and dipole terms. 

Therefore, we may decompose the perturbations Q a B and Kab according to 

c 



QaB = VBqa+£B°Vch a , (21a) 



1 



K AB = 2{V A V B G) TF + -g AB J + 2e (A c V B) V c fc, (21b) 
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where q = q a dx a and h = h a dx a are angular-dependent one-forms on M, and G, J and k are angular-dependent 
functions on M. The advantage of this decomposition relies in the fact that it is covariant. In particular, the 
perturbations arc fully determined by quantities which transform like scalar fields with respect to diffeomorphisms 
on the sphere (i.e. a diffeomorphism ip which leaves M invariant). Here, J = g AB Kab represents the trace of Kab- 
Similarly, we decompose t\b = vs/ + £s Vc<7- In order to define q, h, G, J, k, / and g uniquely, we suppose 
that these quantities lie in the mono-dipole-free space in what follows. 4 As a consequence of our assumptions, the 
transformations p^|) induce the following transformations for the perturbation amplitudes H := H ab dx a dx b , q, h, G, 
J, k: 

H ab h> #ab + V a 4 + V b £ Q , (22a) 
q a i-> q a + £ a + r 2 V a /, (22b) 
fc a i y K+r 2 W a g, (22c) 
G i ^ G + /, (22d) 

J i ^ J+-r a £ a + 2A/, (22e) 
r 

jfe i ^ fc + 3- (22f) 

There are two commonly used methods for dealing with these transformations. The first, which is known as gauge- 
fixing, imposes conditions on the perturbations amplitudes which fix the gauge functions £ a , / and g. A simple way 
of achieving this is to demand k = 0, G = 0, and q a = 0, which simplifies the perturbations considerably. This gauge 
is called the Regge- Wheeler gauge. The second method, called the gauge-invariant approach, does not impose any 
conditions on the perturbation amplitudes. Instead, one constructs linear combinations of the perturbation amplitudes 
which are invariant with respect to the transformations above. For example, we may define the gauge-invariant one- 
form 



h ( : nv) := K - r 2 V a fc, (23) 

>rms trivially, ha i— > ha . 
Pa '■= qa — ^ 2 V a G transforms like p a h- > p a + £ OJ and then setting 



which transforms trivially, h a mv ^ i— > h^™ 1 ^ . The remaining gauge-invariants are obtained by first noticing that 



H% b nV) ■= H ab -V a p b -V bPai (24a) 

j( mv 1 := J--r a p a -2AG. (24b) 
r 

The advantage of the gauge-invariant approach is that no gauge conditions need to be imposed. Therefore, one 
does not need to worry about the physical results obtained depending on a specific coordinate choice. Since the 
field equations are gauge-invariant, it is clear that the linearized field equations can be expressed in terms of such 
gauge-invariant quantities only. 

The variation of the stress-energy tensor, ST, may be decomposed similarly to the variation of the metric: 

ST ab = r abl (25a) 
ST aB = V BfI „+£ B c v cVa , (25b) 

.,2 



ST AB = r z 



2{V A V B a) TF + ^g AB X + 2e (A c V B )Vc/3 



(25c) 



With respect to an infinitesimal coordinate transformation generated by the vector field X we have 8T n- ST + £xT. 



4 For a detailed discussion on monopole and dipole perturbations in the vacuum case, see Refs. [T3. |22|. where it is shown that they 
correspond to stationary modes, i.e. small changes in the mass or small rotations. For the Einstein- Euler case, see Ref. |17| . 
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Explicitly, this gives 



Tab l-> T ao + £{T a b, 

Ha >-> yU Q + T afc C 6 + r 2 PV a /, 

fa l-> + r 2 PV a ff, 

a H> a + P/, 

A i ^ A + -|v a (r 2 P)£ a + 2PA./-, 

/? ^ /3 + P ff , 



(26a) 
(26b) 
(26c) 
(26d) 

(26e) 
(26f) 



where we have used Tab = r 2 PgAB with P := g AB T A B/{^r 2 ). From this, we can construct the following gauge 
invariants: 



(inv 



(inv 
H-a 

(inv 
" a 

(inv 



^(inv 
j^(inv 



■= T ab - P c V c T ab -2T c(a V b)P c , 

■= Va-Ta b p b -r 2 PV a G, 

:= v a -r 2 P^ a k, 

:= a-PG, 

:= A--|v a (r 2 P)p a -2PAG, 

:= (3-Pk. 



(27a) 
(27b) 
(27c) 
(27d) 

(27e) 

(27f) 



For the explicit calculations below, the following trick will be used: as explained above, in the Regge- Wheeler gauge, 
the perturbations simplify considerably since q a = and G = k = 0. Furthermore, in this particular gauge, we 
have ha = h a , H^™' = H a b and j( mtI ) = J, and r^"^ = T ab , fj, a m = jU a etc. Since the field equations are 
gauge- invariant, it is sufficient to perform the calculations in this special gauge; the results in an arbitrary gauge can 
be obtained by simply replacing h a by h a mv \ H ao by H^ b nv \ J by j( mt ') j T(ib by etc. in the final equations. 



B. Computation of the linearized Einstein equations 

As we have just discussed, it is sufficient to consider the following form for the variation of the metric, 

5g a b = H ab , 

Sg a B = £B C ^ch a , 

&9AB = -jr 2 g AB J. 
With respect to this 2 + 2 decomposition the linearized Christoffcl symbols, 

oT M a ^ = - g^ v (d a 8g [iv + d(s5g av - d v 8g af) - 2T a a p5g<j V ) 
= 2 i^aSgpu + VpSgav - ^uSg a p) , 



6T C ab 
6T C a B 



ST C ab 



6T c aB 



AB 



Is- (V. 



Hbd + ^bH a d — Vd-ffa6 ) j 



2 y 



V B H ad + e B C V c (y a h d - V d /i a - 2— h d 



mc ~cd 

01 AB = g 



gABrr a H ad - -gAB^d{r 2 J) + i{ A c V B )^chd 



4 2r^ 

4 (5 c bVa J + S C A V B J - gA B g CD V D j) + gA B e CD V D ( —h 



(28a) 
(28b) 
(28c) 
(28d) 
(28e) 
(28f) 
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Next, we compute the linearized Ricmann curvature tensor, 



Using the expressions (|28|) . the fact that 8R uua p = 8{g lia R a va $) = g^SR va $ + R° ^apSg^a, and the symmetries of 
the curvature tensor, implying 

SRcdab = -^£cd£ab£ e ^ £ 9h 5R e fgh, 



SRcdAB 
SRcDab 



EcdEABE 



xab~CD 



SRa 



bCD, 



e ab e ef 5R cDe f, 



etc., a lengthy but straightforward calculation yields 



SRabed — —-ZabZcd 



V e V f H ef + DH + kH 



SRaBcd 



:£cd { V B 



e ef rV e -H a 



£B C Vc 



-W a {rF) + -e ef r ae h s 
r r 



SRabCD — —-Sab^CD^F, 

SRaBcD = -\^B^ D H ac + \^bd rr b (y a H cb + V c H ab - V b iJ ac ) - ^V a V c (r J) - T -r ac J 



SRaBCD — T;ECD { ~zEB^^ E 



r 2 V a J -2rr h H ab 



rr b £ ab T + Ah a + 2Nh a 



SRABCD = r 2 £ A BECD 



r a r b H ab - ^AJ - lrr a \7 a J -NJ+\j 



(29a) 
(29b) 
(29c) 

(29d) 
(29e) 

(29f) 



where H := g H ab denotes the trace of H ab and where we have defined 



T := r 2 E« b V a ( % 



(30) 



The remaining components of the variation of the curvature tensor are obtained from Eqs. (|29[) using the symmetry 
The linearized Ricci tensor follows from Eqs. (|29[) using the formula 



SRa/3 = 8{g ,iV R^avp) = g^SR^aup - Rfuxvpg^ g VT Sg a 



Explicitly, this yields 



SR„.h = 



j (V a i7 bc + V b H ac - VciTab) - ^2 ^ H -b + kH ab + I g ab (v c V d H cd + OH-~kH 



2r 2 

SRaB = ^B 



1 



V (o [r*V b) J 



\eb C V c 



1 _ 6 - . 2 _. A(r 2 ) Ah a 

~^E a b V b {r l T) + —^ha ' 



SR AB = -^(V A X7 B H) TF +g AB 
£(A C V B) V c (V a / ia ). 



r 
1- 



V>r 6 # a6 ) - 1 rr a \7 a H - ±A(H + J) + ^D(r 2 J 



Finally, we compute the variation of the Einstein tensor using 



(31a) 
(31b) 

(31c) 



SG a/3 = 5R al 3 - -gapg^SR^ + ^g Q pG^Sg^ + (sg a p - ^gapg^Sg^ 
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where here, G := g^G^ denotes the trace of the Einstein tensor. From this, we first find 



SG ab = {5R ab ) tf + \GH% + X -g ah 
1 



\g CD 8R CD + G cd Hl f d + \g cd G cd {H - J) 



SG a B 

SG AB = (SR AB ) TF 



SRaB + ^Ge B c Vch a , 



-gAB 



g cd 8R cd + G cd H l cd - —^^G C d(H - J) 



?cd rrtf 



P.CD, 



Introducing into this the expressions (|31[) for the linearized Ricci tensor and combining the result with the expres- 
sions (|25p for the linearized stress-energy tensor and the background equations ©, the linearized Einstein equations, 
SG nv = kST^l, yield the following set of equations on the two-manifold M, 



e a b V b {r 2 F)~{k + 2)h a = 2 K r 2 (v a - Ph a ), 
V a h a = 2nr 2 fi, 



and 



2rV {a (r c Hg) - rr c V c H a { - ^AH a { + rr {a V b) H - \ V (o (r 2 V b) j) 



T ab -T (a c HH 



tf 



X7 b H tf - — V 



-2V>r fc ^{) + -(A - 2)H + -(A + 2) J + — V Q (r 4 V a J) 
2 ~„, >,„**■. ~„ ~h „tf 1 ^ „ 1 { „ 1 ~ 



H 
1 - 



\v a j 



kt 2 (g ab T ab - T ab Hf f 



b > 



-V a {r»H a {) + V Q VX J b + ~ 2 ^H - — 2 AH _V Q (r^V a J) 

H 

(inv) 



-k(X- pj + PH) 



-Anr 2 a. 



(32a) 
(32b) 



(33a) 
(33b) 

(33c) 

(33d) 
(33e) 



By replacing h a with h a , v a with v a etc. we obtain the corresponding equations in gauge- invariant form, as 
discussed at the end of the previous subsection. For notational simplicity, we omit the superscript ( %nv > in what 
follows. 

We see that the linearized Einstein equations decouple into two groups: the first group comprises Eqs. (|32[) for the 
quantities h ai v a and (3 and describes perturbations with odd parity, sometimes also called axial perturbations. The 
second group comprises Eqs. (|33p for the remaining perturbations amplitudes and describes even-parity perturbations, 
sometimes also called polar perturbations. 



C. The linearized divergence law 



Before we proceed with analyzing the equations, we work out the variation of the divergence law V^T^ = for 
the stress-energy tensor. Using 

<*(V M T^) = d^ST^ + T\ a 8T a fj - r a ^5T» a + T a p 8T% a - T^ a 8T a ^, 

the formula 8T a p = 8{g a ^T IM p) = g afl 8T fl p — g au T^ p8g u)1: and the expressions (|2"5|) we find the following equations: 



1 ^ 



-V a [r> Q - Ph a )] + (A + 2)13 = 0, 



(34) 



in the odd-parity sector and 
1 - 



r V Q (r 2 M a) + (A + 2)a + \{X - PJ + PH) - ]-T ab H ab = 0, 



1 ^ 



(35a) 

V a [r 2 (r ab - H a lT c b )] -^{\-PJ + PH) + ^A/^ - ^T cd V b H cd + ^(T a b - P8 a b )V a J = 0, (35b) 



in the even-parity sector. 

As a consequence of the twice contracted Bianchi identities, V^G^^ = 0, similar equations hold for the linearized 
Einstein tensor. Therefore, the Eqs. (|34I35I) show that the linearized Einstein equations ()32I33|) are not independent 
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from each other. For example, Eq. (|34|) shows that taking the divergence on both sides of Eq. (|32a[) gives an equivalent 
equation than the one obtained by applying the operator A + 2 on both sides of Eq. (|32b|) . Since A + 2 is invcrtiblc 
on the space of mono-dipolc-frce perturbations, we conclude that it is sufficient to impose Eq. (|32a|) in the odd-parity 
sector. Similarly, we see that it is sufficient to impose Eqs. (|33al33bl33c[) in the even-parity sector, the remaining 
Eqs. (|33dl33ep are consequences of the former and the linearized divergence law, Eqs. (|3"5j) above. 



V. VACUUM PERTURBATIONS 

In this section, we analyze the linear perturbations in vacuum. By Birkhoff's theorem (see, for instance, j3lj). the 
background solution must be the Schwarzschild spacetime. The following derivation of the Regge- Wheeler and Zerilli 
equations does not assume particular coordinates on M or S 2 . 



A. Vacuum perturbations with odd parity: The Regge- Wheeler master equation 



Linear perturbations with odd parity in vacuum are described by Eqs. (|32[) with v a = 0, j3 = and P = 0. For the 
following, it is convenient to introduce the coordinate-free notation of differential forms. For a one-form u) = ui a dx a 
on M , we have 



i a b uj b dx a , 



*du = e ab V a uj b 



-V°w a , 



where * and d^ = *d* denote, respectively, the Hodge dual and the co-differential operator on M, respectively. With 
this notation, Eq. (I32a[) reads 



id(r 2 T) - (A + 2)h = 0, 



(36) 

where we recall that T := r 2 e ab V a (r~ 2 h b ) = r 2 *d(r~ 2 h) and where h denotes the one-form h := h a dx a . Applying 
the operator *dr~ 2 on both sides of this equation gives 



12^) 



i(A + 2)^ = 0. 



Setting $ := rJ- and using the background equation Ar = (1 — N)/r = 2m/ r 2 yields the following scalar equation on 
M, 



1 



6m 
r 



$ = 0, 



(37) 



where □ = <$d is the covariant d'Alembertian on (M,g) and m the Schwarzschild mass. Eq. ([37| is the covariant 
form of the Regge- Wheeler equation. Once a solution for $ is known, the metric perturbation h may be reconstructed 
using Eq. ([3"o]) since the operator A + 2 is invertible on the space of mono-dipole-free perturbations. The one-form h 
obtained in this way automatically satisfies Eq. (|32b[) , since Eq. ([3^]) implies that d^h = 0. 



B. Vacuum perturbations with even parity: The Regge- Wheeler and Zerilli master equations 

In the even-parity sector, Eq. (|33ep in vacuum implies that H = 0, and therefore, H a f, = is symmetric and 
trace-free. For the following, it is useful to introduce the one-form 

C := H ab r a dx b . (38) 

As long as iV = r a r a ^ 0, this one- form contains the same information as H ab since 

H a b = (r a C b + r b C a - g a br c C c ) 

if H a b is trace- free. Next, consider Eq. (|33b[) . We may contract this equation once with r a and once with e ah r b . The 
result is, in coordinate-independent notation, 

JtC + ig(dr,dJ) = 0, 
dC - ^dr A dJ = 0, 
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respectively, where we have used the background equation (fob)*-' = 0. The second of these equations motivates the 
definition of the following one-form, 



Z := C -dJ, 
2 ' 



(39) 



in terms of which the above equations read 2d)Z — rDJ = and dZ = 0, respectively. Next, Eq. (|33c[) gives 

4<?(rZ) + r 2 DJ + (A + 2)J = 0. (40) 

Finally, we contract Eq. (|33a|) with r a and use the background equations 2rr a ;, = g a f,(l — TV) and Eq. (j40|) in order to 
eliminate □ J from the resulting equation. This yields Eq. (|41d|) below. 

In conclusion, we obtain the following set of equations governing even-parity linear perturbations of the 
Schwarzschild solution, 



dZ = 0, 

2d t Z + rfij = 0, 

Ad\rZ) + r 2 DJ + (A + 2) J = 0, 



2rg(dr, Z) + 3m J - - (A + 2)(rJ) 



- AZ = 0. 



(41a) 
(41b) 
(41c) 

(41d) 



These equations can be simplified as follows: first, we eliminate DJ from the second and third equation, giving 



1 



2 $Z~2rg{dr,Z) + -(A + 2){rJ) =0 



Using this, we may rewrite Eq. (|41d[) as 



Now there are two ways to proceed. 



d[r 2 <?Z + 3mJ] - AZ = 0. 



(42) 



(43) 



1. We apply the co-differential <v to Eq. (j43|) and use Eq. (|41b[) in order to eliminate DJ from the resulting 
equations. This yields the following equation for 4> := r^Z, 



~ 6m 

-A 

r 



= 0, 



(44) 



which is exactly the covariant form of the Rcgge- Wheeler equation found in the odd-parity sector. 



2. The second way uses the equation dZ = in order to introduce an angular-dependent function (onM such 
that 5 Z = d(. Then, we may integrate Eq. (|43|) and obtain 

r 2 QC + 3mJ- AC = 0. (45) 

Applying the operator (A + 2) on both sides of this equation and using Eq. (|42[) in order to eliminate J yields 



^g(dr,dC)-^A(A + 2)C = 0. 



(46) 



Finally, we define the new scalar ^> by £ = (A + 2 — 6m/r)^ in order to eliminate the first-order derivatives in 
Eq. ([46]) . This yields the following equation, 



1 



6m 



DV+^r A + 2 



j 6m 
(A + 2) z -A+ — 



36m 2 



X „ 2m 
-A - 2 + — 
r 



* = 0. 



(47) 



5 Notice that for the maximal Kruskal extension of the Schwarzschild solution, M is topologically equivalent to R 2 , so this function exists 
according to Poincare's lemma. 
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This is the covariant form of the Zcrilli equation. It involves the inverse of the linear operator A(r) := — A — 2 + 
6m/ r which is invertible on the mono-dipole-free space. Integrating Eq. (|41d[) , one finds the following explicit 
expression for the scalar in terms of the metric perturbations H ab and J: 



* = A" 



A 



6m 
r 



2rr a r b H nh - r r"V„J - -J 



(48) 



This is the covariant generalization of the Zerilli-Moncrief function, see [l5| and references therein. 



Although the Zcrilli equation is more complicated than the Regge- Wheeler equation, its advantage relies in the fact 
that the metric perturbations H ab and J can be reconstructed from the scalar 'J without solving additional differential 
equations. Indeed, from we can compute £ = (A + 2 — 6m/r)$, and knowing £, we obtain the one-form Z = d£ 
and the scalar J from Eq. (|4"2")l . From this, one obtains C = Z + rdJ/2 and then H ab . Using Eq. (|46p in order to 
eliminate □£ the explicit result of these operations, first given in Ref. (32j . is 



J 



H a b 



2 A 



6m 



6m 



1 



2g(dMC)--AC 



2 A + 2 (VaVt^K). 



In contrast to this, the reconstruction of H ab and J from the scalar = 
involved: Integrating Eq. (|43|) using Z = d£, one first obtains 3mJ = A£ 
Eq. ((31]). This leads to the following differential equation for £, 



12mg(dr,dC) - (A + 2) AC 



6m 



r 2 c^Z introduced in method 1 is more 
— d>. from which J can be eliminated in 



4>-{L + 2)4>, 



which could be used to determine the Zerilli potential £ and the function J from 0. However, for m > the operator 
on the left-hand side has a non-vanishing kernel, consisting of superpositions of functions of the form 

( LM (r,#,<p) = e»^ r *Y LM (#,tp), 

where here u>e, m = (L — 1)L(L + 1)(L + 2)/(12m) are the algebraic special frequencies [33j . r* = r + 2m log(r/2m — 1) 
the Regge- Wheeler tortoise coordinate and Y LM denote the standard spherical harmonics on S 2 . 

VI. METRIC PERTURBATIONS COUPLED TO MATTER FIELDS 

Although the Zerilli approach described in the previous section is well-suited for describing vacuum perturbations 
with even parity, it is unclear whether or not it can be generalized to the coupling of matter fields since it is based 
on the integration of Eq. (|41a|) and the introduction of the Zcrilli potential £. Indeed, when matter fields arc present, 
it might not always be possible to replace the symmetric tensor field H a b by a scalar potential £, as is the case in 
vacuum. 6 Here we discuss an alternative approach, which does not rely on the introduction of the Zerilli potential 
and might be more amendable to the coupling of matter fields. In fact, as we show in Sec. IVIII below, it naturally 
leads to a constrained wave system of equations for the case of the linearized Einstein-Euler equations. 

The method we discuss here is based on the observation that in Eq. (|33a[) the gradient of J can be eliminated by 
means of Eq. (|33b[) and the trace of H can be eliminated taking into account Eq. (|33c[) . This leads to the following 



wave-like equation for the trace-free part, H a ^ := of the symmetric tensor field H a 



C ab [H] = 2k 



l_ cd » 2 - 

T ab — T.9 T cd H a b 

2 r z 



Vb)) + 2r 2 V a Vf,a + \2rr( a V b) a + Sr a na 



tf 



(49) 



where the linear differential operator C is defined as 



Cab{H] 



-\7 c \7 c H ab + - ( r c \7 {a H h)c 



b)c 



^r c V c H ab 



Ar A\ . 

2fc + 2 H ab . 

r r A ] 



Important exceptions to this observation include backgrounds which are static, see for example Refs. [8l4l0l| and [l6| for the coupling of 
gravity to linear or non-linear electromagnetic fields. 
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In deriving this equation, we have used the identity 

2V (a (v c # 6)c )]' = (V C V C - 2k) H ab , 



which is valid for sufBcicntly smooth symmetric, trace- free tensor fields H a b- 

The operator C has an interesting symmetry. In order to describe it, we introduce the left-dual of H a b, defined as 7 

(*H) ab := e a c H cb . 

This dual operator maps the space of symmetric, trace-less tensor fields onto itself, and it is invertible since (**H) ab = 
H a b- Now, coming back to the symmetry of £, it can be shown that this operator is invariant with respect to *, that is, 
£ ab [*H\ = (*C)ab[H\ for all sufficiently smooth, symmetric, trace- less tensor fields H ab . Therefore, it is convenient to 
expand H ab in a basis which is invariant with respect to *, since then the equations for the corresponding coefficients 
decouple from each other. 

Such a basis may be constructed from two future-directed null vectors k and 1 with the relative normalization 
k a l a = —2. With respect to a suitable orientation, the metric tensor g ab and the volume form e ab on M have the form 

9ab = — k(Jb)} £ab = +k[J b ], 

and it follows that (*k) a ■= k b e b a = k a , (*l) a '■= lb£ b a = —la- Therefore, k and 1 are eigenvectors of * with eigenvalues 
+1 and — 1, respectively, and we see that k a k b and l a l b form a basis for the symmetric, trace-free tensor fields on M 
which are invariant with respect to the left-dual. As a consequence, we can expand 

H ab = Ck a k b + Dl a l b , (*A)ab = -Ck a k b + Dl a l b . (50) 

The null vectors k and 1 are not unique, since they can be rescaled by a function x on M according to 

kHe*k, ln-e^l, (51) 

corresponding to a local boost with hyperbolic angle x- Following the GHP formalism [HI we say that a quantity / 
has boost-weight q if it transforms like / i-> e qx f under (|5"Tj). For instance, in Eq. (|50]). C and D have boost-weight 
—2 and 2, respectively. 

In order to compute covariant derivatives of H ab we first note that there exists a one-form rj = rj a dx a such that 

Vafcf, = i] a k b , VJ b = -r)Jb, 

due to the fact that k and 1 are null, and due to their relative normalization k a l a = —2 = const. With respect to 
the local boosts, Eq. (fBTj) . the one-form rj transforms like a gauge potential: rj n- rj + dx- Its invariant part, drj, 
determines the Gauss curvature of (M, g): k = —*dr] = —e ab V a rj b . The following "creation" and "annihilation" 
operators 

a+ := fc a (V a - qrja), 

map a quantity of boost-weight q to one of boost- weight q + 1 and q — 1, respectively Furthermore, these operators 
satisfy the commutation relation 

a t-l a q - a q+l a t = 2( fi. (52) 

Now we are ready to compute the covariant derivatives of H ab . Applying V c on both sides of Eq. (|50[) we find, for 
instance, 

V c -ff a b = — i [(aZ 2 C)k c k a k b + (at 2 C)l c kakb + (a 2 D)k c lJ b + (a 2 D)l c l a l b ] . 



One can also consider the right-dual of H ab , defined as (H*) a b '■= H a cC c b- However, it is easy to prove that (H*) a i, = —(*H) a b, so one 
does not obtain a fundamentally new transformation. 
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This immediately yields V a H ab — (a_ 2 C)kb + (a 2 D)l b , and using the commutation relation (|52l) we also find 

(-V C V C + 2k)H ab = {aZ 1 at 2 G)k a k b + (a+aZD)l a l b . 

Using this and expanding — 2r a = r\k a + r^l a , := k a r a , r\ := l a r a , we find the following decomposition of the 
operator C: 



£ab[H\ = 



aI 1 at 2 C + 3-(at 2 C) - — (aI 2 C) + (2— - ) C 



k a hb 



a+aZD + 3— (a^D) - -(a+D) + [ 2— - -_- ] D 



Inlh 



(53) 



Note the symmetry knl which implies C f-> D, a+ -B- a_ 



A. Evolution equations (even-parity sector) 



The perturbation equation (|49| and the decomposition (|53|) of the operator C yield two wave-like equations for the 
quantities C and D, 



aZ^C + 3-(a± 2 C) - — (aI 2 C) + ( 2— - 4 ) C = «T, 



a+a^Z? + 3y (oa -D) - -(4-°) + ( - "a J 15 = kA ' 



(54a) 
(54b) 



where 



T := l a l b 
A := k a k b 



\r a b - \g cd T cd H ab - -^V a (r 2 fi b ) + r 2 V a V h a + 6rr a V b a + 4r a r 6 a 
2 4 v ' 

T, T ab - -7g cd T cd H a b nV„ (r 2 fi b ) + r 2 V a V b a + 6rr a V b a + Ar a r b a 

A 4 r 



(55a) 
(55b) 



When formulating Eqs. ([54]) as a Cauchy problem, one introduces a foliation of M by spacelikc hypersurfaces E t . 
Let u be the future-directed timelike unit normal to and let w be a unitary spacelike vector field orthogonal to 
u, which is, therefore, tangent to E t . We choose the orientation of w such that e ab = — 2u\ a wu. Then, we define the 
null vectors k and 1 as 



k := u + w, 



1 := u — w. 



(56) 



Notice that from this moment on, we lose the boost-invariance leading to the symmetry described in Eq. (|5ip . since 
the foliation singles out a preferred normalization for k and 1. With respect to the 1 + 1 decomposition induced by 
the foliation, the perturbation equations (|54j) can be written explicitly as 



DC 



DD 



^_( 4 Il-2^-4,)c7+(4^-2^-4 M )C' 



r r 

2 (z> - //) - 4 (n 2 - v 2 ) + 8 - v r - ) + 4 ( v- - ft 
AD 



r r 
--ft— 

r r 



4— + 2- - Au )D - 4- + 2 4/x ) D' 

r r I \ r r 



2 (// - v) - 4 U 2 - v 2 ) + 8 ( ft- - v- ) + 4 ( fi- - v- ) + 



r r 

v— 

r r 



8m 



C = K[T + 2{P-~g ab Tab)C] , (57 



D = k [A + 2{P - ~g ab T ab )D] . (57b) 



Here, a dot and a prime refer to the directional derivatives along u and w, respectively, and /i := w a rj a and v := u a r\ a . 
In deriving these equations, we have used the expression k = —s ab V a ilb = fi — v' + fi 2 — v 2 for the Gauss curvature 
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of (M, g), as well as the background equations k = Ar/r — kP and 2Ar/r = 4m/r 3 + K,g ah T a b which follow from 
Eqs. (|6bl6c[) . Finally, m is the Misner- Sharp mass [35| which is defined by r a r a = N = l — 2m/r. Notice the symmetry 

C^D, v^-v, (...)' >-> -{■■■)', T^A, (58) 

which transforms Eq. ()57aj) to Eq. (|57b)l . 

Provided that the source terms Y and A do not contain derivatives of C and D of order higher than one, and 
provided that the quantity J can be eliminated in these terms, Eqs. (|57[) constitute a wave system that, together with 
suitable evolution equations for the matter variables, allows to determine the evolution of the metric perturbations C 
and D from initial data for (C, D) and their first time derivatives on an initial time slice. 

Notice that the wave equations ([57]) decouple from each other for vacuum perturbations of a Schwarzschild black hole. 
In this case, an alternative, coordinate-independent way of fixing the boost freedom in the null vectors k and 1 is by 
specifying the two quantities r\ and rk, which determine the expansion of the two-spheres along in- and outgoing radial 
null geodesies. These two quantities are constrained by the relation rkn = k^ a l b 'r a Tb = —g ab r a rb = —N = 2m/r — 1. 
Note that rk is positive in regions I and III of the Kruskal spacetime and negative in regions II and IV, whereas 
n is negative in regions I and II and positive in the remaining ones. 8 Therefore, the quantities rk and r\ cannot 
be arbitrarily chosen. The connection coefficients ?yk and 771 are obtained from rk and n by solving the background 
equations (|6al6bp which can be rewritten as of (rk) = al 1 (?*i) = 0, a^(rk) = al 1 (?*i) = — 2m/r 2 . 

We choose rk = N and r\ = —1, which is well-behaved in regions I and II (our universe and the black hole region) 
of the Kruskal diagram. In terms of ingoing Eddington-Finkclstein coordinates (v,r), for which g = —Ndv 2 + 2dvdr. 
the corresponding null vectors are 

fca A = 2 7r +N 7T' r TlT = -T-< k - dxa = ~ Ndv + 2dr > Ldx a = -dv, 
ox a ov or ax a or 

and it follows that r/k = 2m/r 2 and 771 = 0. With this choice, Eqs. (|54[) reduce to 

□^-A^ + W-fs-— )A = 



^ , ^k ■ 
r z r r \ r 



with C k = k a V a C and Q = l a V a C. 

A different choice which is valid in regions I and III (our universe and the white hole region) of the Kruskal diagram 
is ^k = >/2) n = —N/y2 which, taking into account changes in the notation and conventions, corresponds to the one 
made in Ref. p9j . where general retarded solutions where computed. With this choice one obtains from Eqs. (|54p their 
Eqs. (2. 21a, 2. 21b) with the identifications C <H- r 2 \l/_2 and D <H- r 2, &2, where ^±2 refer to the complex Weyl scalars 
^2 := R a /3' 1 sk a m^ k 7 m 5 and := Rafi 1 sl a fhPr i m b computed from a Newman-Penrose null tetrad {k, 1, m, m}. 9 
Assuming that the background tetrad is adapted to the spherical symmetry in the sense that k and 1 are orthogonal 
to the metric two-spheres, it is relatively simple to compute the first variation of ^±2 in our perturbation formalism, 
based on the background expressions © and Eq. (|29d[) for the first variation of the mixed components R a BcD of the 
curvature tensor. The result is 

(5*2 = -\m A m B V A V B ( \k a k b H ab + ik a k b V a h t 
r z \2 

where ih := rm, and the expression for 8^-2 follows from this after mapping k 4 1, m 4 m. Therefore, the 
quantities k a k h R a b = 4D and l a l b H a b = 4C determine the even-parity contributions of the variation of the complex 
Weyl scalars ^2 and x &_2, respectively. 



B. Constraint equations (even- parity sector) 

The evolution equations in the previous subsection were obtained from the linearized Einstein equation (|33ap . 
combined with Eqs. (|33b[) and (|33cp . Therefore, to solve the whole set of linearized Einstein equations, one has to 



8 See, for instance, Fig. 6.9 in Ref. [36| for the Kruskal diagram and the definitions of regions I - IV. 

9 In our convention, they are normalized such that g(k, 1) = —2 and g(m, m) = 2. 
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ensure that the Eqs. (|33b[) - (|33ep are also satisfied. First of all, Eq. (|33e[) is automatically satisfied if we impose 
H = —4Kr 2 a. Next, the Bianchi identities and the zero divergence of the stress-energy tensor imply that Eq. (I33d[) 
follows from Eqs. (|33bl33c[) . as discussed in Sec. IIV Cl As we now show, the remaining equations, Eqs. (|33bl33cp yield 
a constraint equation for the evolution system (|54[) plus an equation which allows to determine J from the knowledge 
of (C, D,C,D) and the matter fields on a time slice E t . 

In order to derive the constraint equation, we first project Eq. (|33bj) along the null directions k and 1, obtaining 



a^ 2 C + —a J = nl a (r 2 V a a — /j, a ), 
a 2 D + -aQ.J = nk a (r 2 V a a — /J. a ). 



(59a) 
(59b) 



Applying the operator a~^ 1 to both sides of Eq. (|59a[) and a 1 to both sides of Eq. (|59b[) . and using the commutation 
relation (j52")l with q = 0, one obtains the relation 



aiidi 2 C - a{a 2 D = -2ne ab (V a fi b - 2rr a V b a) =: 2nG 



(60) 



in which J is eliminated. In order to rewrite this equation as a constraint on the spacelike hypersurfaces T, t of M, 
the second time derivatives (C and D) have to be eliminated in the terms alja^C and a~^a 2 D. For this, we note 



a_ 1 a_ 2 C 



! + 2(V W + r/ w ) at 2 C = a_ 1 at 2 C + 2(V W + r/ w )at 2 C, 



where V w := u> a V a and ?; w := w a r] a , and where we have used the definition of the null vectors k and 1 in terms of 
the vector fields u and w, see Eq. (|56[) . in order to write k = 1 + 2w in the definition of o_j. The first term on the 
right-hand side, aI 1 a^t 2 C, can be eliminated using the evolution equation (|54a|) . This yields 



a + n a + 9 C = 



2 V, 



- -(3rj - r k ) 
r 



i± 2 C)-2—(' 



2ijw ) C - 2 



Ar A 



C + kF, 



(61) 



where we have also used a_ 2 C = a_ 2 C — 2 ( V, 



evolution equation (|54b|) . 



2?7w ) C for later convenience. Similarly, one obtains, using the 



a x a 2 D 



2 V, 



'/v 



-(3r k - n) 



(a 2 D) +2— (V w -2r^ 



Ar A \ „ 
£>- 2 D- 



kA. 



(62) 



Therefore, provided that the matter terms T, A and G do not contain derivatives of C or D along u of order higher 
than one, Eq. (|60p with the substitutions for a^aljC an d a^a 2 D given in Eqs. (|61I62[) . yields a constraint for the 
data (C, D, C, D) and the matter fields on a given surface S t . 

Next, we consider Eq. (|33c[) , where we use the background equation ([pa]) , the identity V a V b -ff Q (, = atialaC + 
a7a 2 D, as well as Eq. (|33b[) in order to eliminate the first and second-order derivatives of J. This yields the following 
equation, 



with 



(A + 2) J = 2r 2 (at 1 at 2 C + a -a 2 D) + 4r(r k a+ 2 C + na^D) - 4(r£C + rf£>) - 4t 2 kF, 



F := i.9 ah T Q6 - T ab ff ab + V> a + 4— Ma + r 2 Da + (A — 2)a - 6rr a V a a. 
2 r 



(63) 



(64) 



Provided the derivatives of J can be eliminated in the matter term F, this equation allows one to determine J as 
a function of {C,D) and the matter fields by inverting the operator —(A + 2). In fact, using the substitutions for 
a -i a -2C an( i a-i^D given in Eqs. (|61I62[) . it is even possible to compute J from the data (C,D,C,D) and the 
matter fields on a given time slice St. 



C. Initial data construction (even-parity sector) 



Here, we discuss a practical algorithm for constructing initial data (C,D,C,D) for the wave system ([5"4"|) . We 
assume for simplicity that the matter terms T, A, F and G only depend on (C, D) but not on their derivatives, and 
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that they are completely determined by these quantities and the initial data for the matter fields. Furthermore, we 
assume that G does not depend on C nor on D. As we will see in the next section, these assumptions are satisfied in 
the perfect fluid case, for instance. 

Our construction specifies J and its time derivative, J := u a V a J, on the initial surface, as well as suitable initial 
data for the matter fields. As a consequence, the full gradient V a J is known on So, and we can determine the 
quantities <iL 2 C and a 2 D from Eqs. ()59a|59b[) . Next, we combine Eqs. (|60I61I62I63|) in order to eliminate ol 1 aI 2 C 
and a^a^D, and obtain the two equations 




\C + r 2 D) 



2 C) + ^(a 2 D) + J^(A + 2) J + «(T - F - G), 



(65a) 



{rlC + rfD 



(a± 2 C) + —{A + 2)J- K {A-F + G). (65b) 

T \T 



The first three terms on the right-hand sides of these equations are known, since a^i 2 C, a 2 D and J are known on 
the initial surface So- Furthermore, according to our assumptions, T, F and G depend on (C,D) but not on their 
derivatives. Therefore, Eqs. (|65p constitute a system of ordinary differential equations for C and D. Once this system 
has been solved, the initial data for (C, D) is obtained by setting 



C 



C - 2rjkC + a± 2 C, D := D' + 2^D + a+D 



(66) 



This provides initial data (C, D, C, D) on the initial surface Sq- The evolution equations (f5"T|) . together with suitable 
evolution equations for the matter fields, determine (C, D) on an arbitrary time slice E t . The metric field J is obtained 
from this by solving Eq. (|63l) . This determines the metric perturbations in the even-parity sector. 



D. 



The odd-parity sector 



Odd-parity gravitational perturbations are much simpler to describe than even-parity ones, and gauge-invariant 
master equations for them have been already obtained in the literature, see for instance Refs. (l2l. |37(. Since our 
approach does not require the decomposition into spherical harmonics, and also for completeness, we briefly describe 
the odd-parity sector here. 

In this case, the metric perturbations are determined by the one-form h on M whose dynamics is governed by 
Eqs. (|32|) . As discussed in Sec. IIV CI it is sufficient to consider the equations of motion for the matter fields and 
Eq. (|32a|) . which, in terms of the notation of differential forms, reads 



id(r 2 T) - (A + 2)h = 2w 2 Y, 



(67) 



with T = r 2 *d(r 2 h) and T := [v a — Ph a )dx a . Now we can proceed exactly as in Sec. IV Al to derive a master 
equation for the scalar field $ := r? ' . Applying the operator r 4 *dr~ 2 on both sides of Eq. (|57)l we first obtain the 



following inhomogeneous wave-type equation for r T: 



- (A + 2){r 2 F) = 2«r 4 ?(iT, 



Introducing T = $/r and using the background equation (|6b[) . we then obtain 



A 



6m 



9 J-ab 



$ = 2KrJdT, 



(68) 



(69) 



which generalizes the Regge- Wheeler equation (|37|) to the case of odd-parity linear metric perturbations coupled to 
matter fields. Provided the one-form Y does not depend on h, this yields a wave equation for the gauge-invariant 
scalar $ on M. Once $ is known, h can be reconstructed from Eq. (157]) on the mono-dipole-free space. In case Y 
depends on h, the method just described also works provided Eq. (1571) can be solved for h. 
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VII. FLUID PERTURBATIONS 

Here, we apply the ideas described in the previous section to linear perturbations of spherical, self-gravitating 
fluid configurations. We start in Sec. IVH Al with a short review of the equations of motion for a relativistic perfect 
fluid in thermodynamic equilibrium and discuss potential flows. Next, in Sec. IVIIB1 we linearize the flow, assuming 
that the background is vorticity-free and isentropic, and show that the fluid perturbations can be split into two 
parts, the first part describing the perturbations of the vorticity and the entropy and the second part describing a 
linearized potential flow. The resulting equations are specialized to the case of a spherically symmetric background 
in Sec. IVII Cl In Sec. IVII Dl we linearize the stress-energy tensor and the propagation and constraint equations in the 
odd- and even-parity sectors are derived in Sec. IVII El and Sec. IVII Fl respectively. 

Previous work regarding the linearization of relativistic perfect fluids include Rcf. [38] , where the stability of steady 
state accretion flows is analyzed and on which most of the ideas presented in Sec. IVII Al and Sec. IVII Bl are based, 
Ref. [1^] which analyze the stability of static stars and Ref. [I?} which derive a system of propagation equations based 
on the covariant gauge-invariant approach. Besides the avoidance of the spherical harmonic decomposition, the results 
presented below differ from the ones obtained in Ref. [TtJ in the strategy for deriving the evolution equations. While 
in (l7| projections of the linearized stress-energy tensor which eliminate the matter variables are considered, here we 
exploit the fact that a spherically symmetric isentropic flow gives rise to a closed one-form v, which allows to naturally 
split the fluid perturbations into a linearized potential flow and a complementary part describing the propagation of 
the linearized vorticity and entropy. Although the results in [l7| lead to a free (unconstrained) evolution system and do 
not require the background to be isentropic, we believe that our approach allows for a clearer physical interpretation 
of the fluid perturbations. 

A. Basic fluid equations and potential flows 

We consider a perfect, relativistic fluid with isotropic pressure. It is described by the energy density e, the pressure 
p, and the particle density n, as measured by an observer co-moving with the fluid elements which have four-velocity 
u (u is normalized such that g^u^u" = — I). The equations of motion for the fluid are 

V M J" = 0, V„T^ = 0, (70) 

where J M = nu^ is the particle current density and T^ v = eu^u v + pj^, 7^ := g^ v + u fJ, u L ' , denotes the stress-energy 
tensor. We require local thermodynamic equilibrium, that is, each fluid element is in thermodynamic equilibrium. 
Therefore, we assume the existence of an equation of state e = e(s, n) which satisfies the first law of thermodynamics, 

d(i)=T*-p*(i), 

with T the temperature and s the entropy per particle. For the following, it is convenient to replace the energy density 
e by the enthalpy h := (e +p)/n per particle, which satisfies dh = Tds + dp/n. With this notation and under these 



assumptions the fluid equations (|70[) yield (assuming T > 0), 

V u s = (entropy conservation along the fluid lines), (71a) 

V u n + On = (particle conservation law), (71b) 

ha^ + D^h — TD^s = (relativistic Eulcr equations), (71c) 



where 9 := V ^u^ and := V u w^ denote the expansion and acceleration of the fluid, respectively, and where 
Dfih := 7 / j 1/ V„/i denotes the differential of h projected onto the space orthogonal to u. Equations (fTTj) . together with 
the equation of state form a closed evolution system for the quantities s, n and u. Eq. (|71a[) implies, in particular, 
that there is no heat transfer between the different fluid elements since Tds — along the fluid lines. (However, the 
entropy s may vary from one fluid trajectory to another.) 

In terms of the one-form v := hu^dx^ and the vorticity f2 M „ := 7 M Q 7„ .^V '[ a u p] °f the fluid, the Euler equations (|71c|) 
are equivalent to 

F^ v := V^tv - VuV^ = 2 {hVt^y - Tu^D^s) . (72) 

Therefore, the one- form v is closed if and only if the flow is irrotational = 0) and isentropic (s = const). In this 
case, there exists (at least locally) a potential ip such that v = dip. Then, the particle conservation law, Eq. (|71bp . 
yields 



v M (£ V ^) = 0. 



(73) 
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Here, n is regarded as a function of h obtained by inverting the relation h = h(n) and h is obtained from the condition 



n — h 2 u^u^ = —h 2 , which yields 



h = y/-V*il> ■ V„V- (74) 



For the case of a stiff equation of state, for which h is proportional to n, Eq. (|T3[) reduces to the standard wave equation 
O g tp = on [M, g) for the potential ip. However, in general Eq. f|T3f) is a nonlinear wave equation whose characteristics 
depend on the gradient of ip. Eq. ((75)) and the associated stress-energy tensor T MU = nhr 1 (V ,,J "ip)(W' 'ip) +p(h)g fJ,L ' can 
also be obtained from the simple action functional 

S^, E ] := J p(h)^\f\d 4 x, (75) 

where h is given by Eq. (|74j) and p(/i) obeys the first law dp = nd/i. 

The linearization of Eq. ([75)1 yields the remarkably elegant equation [38| 

□«(5V = (76) 
for the perturbed potential Sip, where □© is the wave operator belonging to the sound metric 25 defined by 

h v s oe olog(n) 

w s denoting the sound speed (in units for which c = l). 10 For the following, we assume that the sound speed satisfies 
< v s < 1. By employing an ortho normal tetrad {eo,ei, 62,63} with respect to g which is adapted to the fluid flow, 
i.e. such that eo = u, it is not difficult to see that © is Lorentzian, and that its cone (the sound cone) lies inside or 
(if v s — 1) on the light cone of g. Notice also that u is timelike with respect to both g and ©. 

B. The linearized fluid equations 

Linearizing Eq. (j72|) about an iscntropic, irrotational fluid configuration one obtains 

V M (*u„) - V„(^) = 2 (hSQ^ - Tu {ll D v] 8s) . (78) 

For the following, it is convenient to decompose Sv into two terms, the first one being the differential of a function \1/ 
and the second one a one-form W which is orthogonal to u: 

Sv = d^ + W, u^W^ = 0. (79) 

This decomposition always exists. Indeed, if (Sv is given, we may integrate the equation £ u ^f = u^Sv^ along the flow 
lines to obtain "I 1 , and set W := Sv — d^>. By construction, u M W^ = and Sv = dty + W. The quantities \& and W 
are unique up to transformations of the form 

* m. * + /, Wy. h+ Wy. - DJ, (80) 

where / is a function satisfying £ u f = 0. This freedom can be fixed by restricting the initial data to 1 J f = 0, or by 
imposing suitable initial conditions on W. 

Introducing the decomposition (|79|) into Eq. ([75]) and projecting the result orthogonal to u gives 

2h5^ u = D ll W v -D v W^ (81) 

which shows that the one-form W describes the vorticity part of the perturbations. The remaining information in 
Eq. (|78|) is obtained by contracting with u^. Using it M W^ = 0, this yields 

£ U W. = TDJSs), (82) 



10 The following expressions for the inverse sound metric and its determinant are useful: 



<6» v = -vs 
n 



vW\ = Vlsl 
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which, together with the linearization of Eq. (|71a[) , 

£ u (6s) = 0, (83) 

describes the evolution of the variation of the entropy and the one-form W along the flow. As a consequence of 
Eqs. (|81|82J) one also finds 

£ u {h6% v ) = hDfr D v] 5s, (84) 

which shows that in general, a spatial gradient of 6s generates vorticity along the flow. Eqs. (|82|) and (|83| describe the 
evolution of the variations of the entropy and the one-form W which describes the variation of the vorticity through 
the relation (|81j). These equations form a system of advection equations along the fluid lines, and since they do not 
depend on \ff nor on the variation of the metric fields, they decouple from the remaining perturbation equations and 
can be solved separately 

It remains to linearize Eq. (|71bp . which is equivalent to the continuity equation V tl (nh~ 1 g^ v v w ) = 0. The variation 
can be written in the form 



v„r = 0, 



where the vector field 3 M is given by 



(85) 



In order to evaluate the variation of n/h, we regard n = n(h, s) as a function of the enthalpy h and the entropy s per 
particle, which is obtained by formally inverting the relation h = h(n, s). Using the definition of the sound speed v s 
and the first law, we find the thermodynamic relations 



dp 




d log h 


dh 


de 


s 


dlogn 





n os 



Using the identities 



dp 


ds 


dv 


ds~ 


v dv 


P Q P 



dv 




dv 


dT 


ds 


p 


df 


P ds 



v := -, 
n 



and the definitions of the heat capacity per particle at constant pressure c p , the coefficient of thermal expansion a p 
and the compressibility at constant entropy K s , defined as (see, for instance, Ref. |40| ) 



T to 
dT 



1 dv 
v df 



1 dv 

v dp 



the term involving the partial derivative of the pressure can be rewritten as 



Gathering the results, we obtain 



dn 
dh 



n 1 

h vl 1 



from which we can finally compute 



dn 
~d~s 



op 



= T- 



p 



CpK, s 



n 1 dh 

h v? ds 



n T 
h vl 



A, 



n 

h 1 



a/3, 



AT 



A := 1 



1 a r , 



Tit Cp 



(86) 



(87) 



Next, we compute the variation of the identity h = —g a ^v a v^ which follows from the definition of w M = hu^ and the 
normalization of u. Taking into account the decomposition (|79l) this gives 



2h6h 



-2v a 6v a + v a v^6g a p = -2u Q V Q * + v a v f3 6g a0 . 



(88) 
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Combining Eqs. (185I87I88[) and recalling the definition of the sound metric, Eq. (l77l) . we finally arrive at the following 
expression: 

W (V,* + W v - Sg^v 13 ) + \v^ a0 8g afi - —ATds . 
v '2 nv s 

Therefore, the linearized continuity equation yields the inhomogeneous wave equation 

□<S# = div e (X), (89) 

where div© refers to the divergence of a vector field with respect to the Levi-Civita connection associated to the sound 
metric, and X is the vector field defined as 

X" = <&» V {W V - 6g v pvP) + \v^ af3 5g a p - —ATSs. 

2 nv s 

Notice that X vanishes when W, Ss and Sg are zero, and in this case Eq. (|59")l reduces to the homogeneous wave 
equation (|76p describing a linearized potential flow. In the general case, the one-form W and the variation of the 
entropy, Ss, are obtained by integrating the advection equations (|82I83[) . and therefore, they can be considered to be 
known quantities in the expression for the vector field X. However, the variation of the metric fields, (Sg, are coupled 
to the potential \f r through the lower order terms on the right-hand side of Eq. (|59")l and through the linearized Einstein 
field equations. 




C. Fluid perturbations on a spherically symmetric background 

So far, we have only assumed the background to be vorticity-free and isentropic. From now on, we assume in 
addition the background to be spherically symmetric. In particular, this implies that the four-velocity vector u has 
no angular components; hence it can be regarded as a vector field on the two-dimensional manifold M and our 
zero vorticity assumption follows from the spherical symmetry. Defining [I?} the additional vector field w whose 



components are w 



-e a bU , we obtain an orthonormal basis of vector fields satisfying g a b 



-u a u b + w a Wb and 



i a b = —2u[ a Wb]. The connection coefficients with respect to this basis are entirely determined by the two quantities 
/i := V a u a and v := V a w a . 

Since the one-form W is orthogonal to u, we can expand it according to 



W a = Uw a , 



W A - V A V + i A B V B u, 



where U, V, and uj are angular-dependent scalars on M, and where we have used the decomposition described 
in Eq. (JTTJ. Performing the calculations in the Regge- Wheeler gauge, for which 8g a b = H a b, Sg a B = £b vch a , 
SgAB = r 2 gABJ/2, we first obtain from Eqs. (|85I87[) and 



5h 
~h 

r = 



l£ u m + ^H ab u a u b , 



ii 

Ti 

l n r 



+ Uwb 



H bc v c )+-(a cd H cd 



J)v a 



A 



^ r 2 h 



V A {^ + V) + e AB ^J b {lo - h a v a ) 



where we have defined the sonic two-metric 



~ ab 



Therefore, the linearized continuity equation gives 



r h 



Uw a - Q ab H bc v c 



^v a (a bc H bc 



J)-u a ^TSs 



(90) 



-AV. (91) 



The right-hand side of this equation could be further simplified by using the background equation V a (r 2 nu a ) = 
V a {r 2 nh~ 1 v a ) = 0. For a spherically symmetric background, Eqs. (|52")l and (|83p imply the ordinary differential 
equation 



£ u u> = 0, 



(92) 
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in the odd-parity sector, and in the even-parity sector the system of differential equations 

(£ u + n)U = T(5s)', £ U V = TSs, £ u (Ss) = 0, (93) 

with (6s)' = w a \7 a (Ss). 

In the next three subsections we combine the above results with the equations describing gravitational perturbations 
derived in Sec. IVI[ obtaining effective equations describing the linear fluctuations of a spherically symmetric self- 
gravitating fluid. 



D. The linearized stress-energy tensor 

The variation of the stress-energy tensor = nh v^v^ +pg^ v yields the general expression 



ST„ 



9^ 



1 



1 u,,u,, 



1 



A 



TSs + 2nu (lx [V v) * + W v) ] + pSg„ u , 



where we have used Eqs. (|87I88I79[) and the hrst law of thermodynamics dp = ndh—nTds. For a spherically symmetric 
background, the corresponding quantities r a b, fx a , v a , a, A and /3 defined in Eqs. ([23)1 are 



v a = PK + nu a u, /3 = 0, 



(94) 



in the odd parity sector, and 



Tab = n 



fjab + I "4 - 1 I U a U b 



]-u c v d H cd - u c V c y 



A 

fjab + -^U a U b 

vi 



TSi 



2nu 



(a 



V b) y + Uw b} 



pHab 



Ha = nu a (V + V), 
a = 0, 



A = n 
in the even parity sector. 



u a v b H ab - 2(u a V a * + TSs] 



■pj, 



(95a) 

(95b) 
(95c) 

(95d) 



E. The linearized Einstein-Euler equations: odd-parity sector 



The odd-parity perturbations of a spherically symmetric fluid configuration are described by the linearized Euler 
equation (|92[) and the linearized Einstein equations (f32j) . with the coefficients v a and j3 given in Eq. (|94|) . which yield 



id(r 2 J-) - (A + 2)h 



2Kr 2 nuuj : 
0, 



(96a) 
(96b) 



where J- = ?* 2 Sd(r _2 h) and u := u a dx a is the one-form on M corresponding to the four-velocity of the background 
flow. We have already noticed in Sec. IIV CI that Eq. (|96b|) follows from Eq. (I96a|) and the divergence law for the 
stress-energy tensor. In the present case, this can also be verified explicitly by applying the co-differential on both 
sides of Eq. (|96a[) , using £ u uj = and the background equation <r (r 2 nu) = 0. 

A wave equation for the scalar field $ := rF follows from the method described in Sec. IVI Dl For this, we use the 
background equation dv = and first note that 



idT = *d(nuuj) = id ( -wv) = * 



h ( — uj 



where the prime refers to the directional derivative along w. With this observation, Eq. yields (cf. Appendix B 
in Ref. (H): 



- (A + 2){r z F) = -2 K r 4 /i -w 
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Introducing T = $/r and using the background equation f|6b[) . wc can rewrite this equation in the form 



A 6m k . 
72- — +2(^P) 



$ = -2nrh ( -w 



(97) 



which generalizes the Regge- Wheeler equation (|37[) to the case of odd-parity linear fluctuations of self-gravitating 
spherical fluid configurations. Once $ is known, the metric perturbation h can be reconstructed from Eq. (|96a[) on 
the mono-dipolc-frce space. 



F. The linearized Einstein-Euler equations: even-parity sector 



The even-parity perturbations are described by the linearized fluid equations (|9 1 193[) and the linearized Einstein 
field equations (|33|) . where the coefficients T a \, and fj, a , a and A are given in Eqs. (|95|) above. First of all, we recall that 
the Eqs. ()93|) are simple advection equations along the flow lines, which are decoupled from the remaining perturbation 
equations. Therefore, specifying initial data for ([/, V, Ss) on an initial hypersurface So completely determines the 
evolution of these quantities. Next, as discussed in Sec. IVI A[ the linearized Einstein equations give rise to a system 
of two wave equations for the metric fields (C, D), see Eqs. (|57j) . where the source terms T and A are determined by 
the linearized stress-energy tensor. Here, the vector fields u and w which were introduced in Sec. IVI Al to determine 
the null vectors k := u + w and 1 := u w are oriented such that u coincides with the four- velocity of the background 
fluid flow. Using the definitions in Eqs. ([55|) and the expressions (f95|) . we find 



1 

-2-1 



-- 1 

vi 



D)-i$> 



2 v 2 s 
2 vi 



ihC 



nhD 



— Vi (r 2 n) - (fi - v)n 



(r 2 n) 



(// + v)n 



(* + V) + n fViV + U 

(* + V) + n ( v k y - u 



Notice that Eq. (|57a[) is transformed into Eq. (|57b[) one under the symmetry ([55)1 



where A is defined in Eq. 
augmented by U n- —U. 

Together with the evolution equation for 'J, Eq. ([9"Tj) . we obtain a closed wave system for the three quantities 
(C, D, \E'). In order to write down the evolution equation for VP, we notice that J can be eliminated from the right- 
hand side of Eq. (f9"Tj) by virtue of Eq. (|33b[) and the satisfaction of the background equation V a (r 2 nu a ) = 0. This 
yields the equation 



7 a (fl^Vb*) 



V a log( — 



V k (r 2 n) 



g a6 V 6 * - -^(tf + V) - 2«n/i(* + V) 



2A 



__i)(i7 + ij)__r«y a 



-ijt + v) 



C + h 



Vi(r 2 n) 



- (A* - v) 



D + U' + vU 



log 



£/, (98) 



describing the evolution of the fluid field <£. Here, g ab refers to the two-dimensional sonic metric on M, see Eq. (|90l) . 
The evolution system (|57I98[) for (C, D, vf) is subject to the restriction in Eq. where the quantity G is 

G= [n(* + V)]' + m>(V + V), 

for the fluid case. Notice that G does not depend on J; hence, Eq. (|60|) gives a relation between (C, D) and the matter 
fields (^>, V). The quantity J can be reconstructed from the variables (C, D,^,U, V, 5s) using Eq. (f6"3")) where 
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'■(C + D) 
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^j-V u (r 4 n) + /m 



Initial data for the evolution system (|57I98[) can be constructed by specifying (J, J, (7, V, <5s) on an initial 
hypersurface E , and by following the algorithm described in Sec. IVI Cl Namely, we first set 



a+ 2 C 



a 2 D 



a J + fcn(* + V), 
a^J + Kn(* + V). 
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Next, we solve the ordinary differential equations (|65j) for (C,D), where the matter source terms are given by the 
expressions r - F -G = X -Y + nhC, A - F + G = X + Y + nhD with 

+ n(l-^j TSs - 2n + (A (* + V), Y := n{^' + 2V 1 - U) + ^(rn)'(tf + V). 

Finally, (C, D) are determined by Eq. (|66|) . This provides initial data for the wave system (|57I98[) . At this point, we 
recall the freedom described in Eq. (|80| that allows to set the initial value for ^ to zero. Therefore, our initial data 
contains six degrees of freedom, in accordance with the results in fl7| . 



X := 



h 



(C + D) - * 



VIII. SUMMARY AND CONCLUSIONS 



In this paper we have presented a gauge-invariant perturbation formalism for spherically symmetric background 
configurations in general relativity, without assuming the background to be static or vacuum. The new feature of this 
formalism is that it combines the covariant, gauge- invariant method in fl2j with the quasilocal approach in [22j, so 
that it does not require a decomposition of the perturbations into tensor spherical harmonics. 

As a pedagogical first step, we have derived effective wave equations on the two-dimensional orbit manifold 
A/ = M/SO(3) orthogonal to the two-spheres, describing the propagation of scalar and electromagnetic test fields on 
an arbitrary spherically symmetric background. Then, we have further developed our formalism to describe metric 
fluctuations. Our approach is based on the 2 + 2 form of the background geometry and the construction of a full 
set of gauge-invariant amplitudes. These amplitudes are angular-dependent tensor fields on M, which behave as 
scalars under rotations of the metric two-spheres. We have expressed the linearized Einstein equations in terms of 
these gauge-invariant quantities which decouple into two groups, describing perturbations with odd and even parity, 
respectively. 

Next, we have applied our formalism to the vacuum case and have derived the covariant generalizations of the Regge- 
Wheeler and Zerilli master equations, which describe the propagation of arbitrary linearized gravitational waves on 
a Schwarzschild black hole. In particular, we have shown that the Regge- Wheeler equation can easily be obtained in 
both, the odd- and the cuen-parity sectors, although in the even-parity sector it does not represent a master equation 
in the strict sense, since in this case an additional differential equation needs to be solved in order to reconstruct the 
metric perturbations. 

For a Schwarzschild background, the master equations describing the propagation of scalar, electromagnetic, and 
linearized gravitational perturbations have the form of an effective wave equation on M, 

□* + V(r)$ = S, V(r) := \ (-A + fi(r)) , (99) 

where □ is the covariant d'Alambcrtian on the radial part (M, g) of the Schwarzschild-Kruskal manifold, A denotes 
the Laplacian on the sphere S 2 and the operator B(r) is defined as 

for a Klein- Gordon field 

for odd- and even-parity electromagnetic fields 
for odd-parity gravitational perturbations 
) A(r)~ 2 , for even-parity gravitational perturbations 

with A(r) := — A— 2+6m/r. Here, the source term S is zero for the Klein-Gordon field, and in the electromagnetic case 
it is determined by the four-current charge density, see Eqs. (|14cll6p . In the gravitational case, S is zero for vacuum 
perturbations. For gravitational perturbations generated by an infinitesimal stress-energy tensor, the expressions for 
S are given in [TH]. Recall that in the electromagnetic and gravitational cases, $ is restricted to the monopole-free and 
the mono-dipolc-free spaces, respectively. This means that the expansion of $ in spherical harmonics only contains 
terms with angular momentum number I > L with L = in the scalar, L = 1 in the electromagnetic and L = 2 in the 
gravitational cases. Based on the energy-type estimates described in Q, one can prove that the solutions of Eq. (f99|) 
belonging to sufficiently regular initial data on a spacelike slice remain uniformly bounded outside the black hole. 

Next, we considered metric perturbations coupled to arbitrary matter fields. While in the odd-parity sector the 
generalization to matter is rather straightforward in most cases, in the even-parity sector it is unclear whether or 
not the Zerilli approach works. Recall that in vacuum, Zerilli's method in the covariant approach is based on the 
observation that one of the linearized Einstein equations implies the zero exterior derivative of the Zerilli one- form Z, 
see Eq. (|41a[) , which means that Z can be written as the differential of a scalar potential. However, with the exception 
of particular matter models, obtaining a closed one-form does not seem possible in the general case. Therefore, we 



B(r) := { °'6m 

r ' 

6m _ 24m / ^ _ 3m\ j^/j-)- 1 
r r \ r J ^ ' 



12m- 1 

—?r- 
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proposed an alternative approach which naturally leads to a constrained wave system for two gauge-invariant metric 
perturbation amplitudes. In vacuum, this wave system decouples into two wave-like equations which are related to 
the Teukolsky equations [28l. |29| for the Weyl scalars and ^4. In the presence of matter fields, the two equations 
are coupled to each other through the linearized stress-energy tensor, and, together with suitable evolution equations 
for the matter fields, describe the propagation of the linearized gravitational and matter field perturbations. 

As an application, we considered linear perturbations of a spherically symmetric perfect fluid configuration. As- 
suming local thermodynamic equilibrium, and assuming that the background is isentropic, we decomposed the fluid 
perturbations into two parts. One part describes a linearized potential flow and is entirely determined by a scalar 
quantity The complementary part is described by a vector field W which is orthogonal to the four- velocity of the 
background flow, and generates a non-trivial vorticity field at the linearized level. The evolution of W couples only 
to the evolution of the perturbed entropy, and the corresponding equations are simple advection equations along the 
background flow. Then, we combined these equations with those describing gravitational perturbations and obtained 
effective equations on M describing linear fluctuations of a spherically symmetric self-gravitating fluid. In the odd- 
parity sector we obtained an inhomogeneous master equation generalizing the Regge- Wheeler equation to the fluid 
case where the source term is determined by the fluid field W. In the even-parity sector, we obtained a constrained, 
coupled wave system for two gauge- invariant metric perturbation amplitudes and the fluid potential \I r . This evolution 
system may be solved using numerical methods, regarded as a Cauchy problem with initial data constructed according 
to the algorithm given in Sec. IVI CI 

The covariant propagation equations found in this article should have many applications in astrophysics. For 
instance, we can mention the stability analysis of radial accretion flows under linear perturbations. Another important 
application is the stability analysis of Cauchy horizons found in some spherically symmetric models, like the Rcissncr- 
Nordstrbm black hole or a spherical dust cloud undergoing complete gravitational collapse. The latter is particularly 
interesting, since it may lead to the formation of globally naked singularities which arc stable with respect to spherical 
perturbations [ill |42|. 

The fact that our formalism eliminates the need of decomposing the perturbations into spherical tensor harmonics 
represents a simplification with respect to previous approaches, which could be useful for a generalization to second- 
or higher-order perturbations. For example, the equations governing the second-order perturbations contain source 
terms which depend quadratically on the first-order ones, and when expanded into spherical harmonics this leads 
to the computation of Clebsch-Gordan coefficients. In our approach, this computation is unnecessary. Instead, 
the source terms need to be decomposed according to Eqs. (|A2IA3[) which can be performed by solving the elliptic 
equations (|A4IA5|) on the two-spheres. 

Another potential advantage of the results derived in this article is related to the constrained wave system describing 
the gravitational perturbations in the even-parity sector. Since one of the scalars involved in this system is directly 
related to the Weyl scalars ^4, it should be simple to compute the gravitational radiation emitted by an isolated 
system. 

Finally, it should also be possible to generalize our formalism to higher-dimensional, 50(<7+l)-symmetric spacetimes 
of the form M = M x S q , where the dimension of the invariant spheres, q = 2, 3, . . ., is general. As long as q = 2, 
we can base the construction of gauge- invariant tensor fields on M on the same decomposition as in Eqs. (|A2IA3p . 
although the effective equations obtained on M are now more complicated if the dimension of M is higher than two. 
When q > 2, one has to consider the decompositions (IA6IA7I) instead, and in this case, the equations decouple into 
scalar, vector, and tensor perturbations. See Refs. [43l.l44l] forapplications to the stability of higher-dimensional static 
black holes. 



Acknowledgments 

It is a pleasure to thank T. Zannias for insightful discussions and comments and D. Brizuela for reading the 
manuscript. This work was supported in part by Grants CIC 4.19 to Univcrsidad Michoacana and CONACyT 46521 
and 238758. 



Appendix A: Covariant, orthonormal decomposition of one-forms and trace-free symmetric tensor fields on 

compact two-surfaces 

Let (S, g) be a compact, two-dimensional Riemannian manifold without boundary which is orientable and has 
positive Gauss curvature (By the Gauss-Bonnet theorem, this implies that S is topologically equivalent to the two- 
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sphere S 2 .) Let e be the associated volume form on (S,g) and introduce the following scalar products: 

<f,9>:= J f-g, <u>, V >:= J g AB Lu A - VB , < X, Y >:= J g AB g CD X AC ■ Y BD , (Al) 
s s s 

for smooth functions /, g, smooth one-forms u>, rj and smooth tensor fields X, Y on S. In this appendix we show: 

Proposition 1 Let u) and X denote a C°° one-form and symmetric, trace-less tensor field, respectively, on (S,g). 
Then, there exist C°° functions f , g, h, k on S such that 



UA = VAf + e A B V B g, (A2) 
HA C V B} ' 



Xab = {V A V B h) TF + e {A c V B) V c k, (A3) 



where Y AB := Y AB — \gABg CD YcD refers to the trace-free part of the tensor field Y on S. 

Furthermore, the decompositions I(A21A3\) are unique and orthogonal with respect to the scalar products defined in 
\A1\) . The functions f and g are unique up to an additive constant and the constants h and k unique up to the 
addition of an element in the kernel of the operator (VaVb) tf . 

Remarks: 

1. The functions (/, g) in Eq. (| A2|) can be found by solving the elliptic equations 

A/ = V a uj a , Ag = -e AB V A u; B , (A4) 
and the functions (h, k) in Eq. (|A3[) by solving 

A(A + 2k)h = 2V A V B X ABl A(A + 2k)k = -2V A V B e A c X BCl (A5) 
with A = \7 A \7 A the Laplacian and k the Gauss curvature with respect to (M,g). 

2. For the case of the unit two-sphere, S = S 2 , the kernel of (VaVb) tj? consists of the functions of the form 
/ = fooY 00 + fnY 10 + f n Y u + fnY 1 - 1 with foo, f w g E, f n g C, and Y im denoting the standard spherical 
harmonics. In the context of the perturbation theory described in this article, we call the constant part /oo^ 00 
°f / = 2 ftmX im the monopole part, the part of / which is equal to fwY 10 + fnY 11 + fnY 1 " 1 the dipole part 
and the remaining part / = X^>2 fzmY im the mono- dipole- free part. In this case, it is not difficult to prove the 
proposition by decomposing the elliptic equations (|A4[) and (|A5[) into spherical harmonics. 

3. The decomposition (|A2I) also follows from the Hodge decomposition theorem (see, for instance Proposition V.8.2 
in Ref. |45j ) and the fact that all harmonic one-form vanish on (S,g) with our assumptions, see Lemma [1] below. 

Proof of proposition [TJ We split the proof in two parts. In the first part we show there exist C°° functions / and 
h and a divergence-free C°° one-forms rj and a C°° symmetric, trace-free tensor field Y satisfying \7 A S7 B Y AB = 
such that 

LO A =V A f + T] A , (A6) 

X AB = (V A V B h) TF + Y AB . (A7) 

These decompositions are orthogonal with respect to the scalar products defined in Eq. (|A1[) : hence, they are unique 
(if they exist). The first equation implies that / must be such that the one- form u) — V/ is divergence- free, which is 
equivalent to the statement 

< oj - V/, Vu >= (A8) 

for all C°° functions u on S. Similarly, the function h must satisfy the problem 

< X - (VV/i) Ti? , VV« >= (A9) 

for all C°° functions v on S. In order to show that the problems (|A8IA9p possess a solution, we use standard tools 
from the theory of Fredholm operators, see, for instance, Ref. (46[, or Theorem 4 in appendix A of [13] for a summary 
of the relevant results. 
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Denote by H k (S) the Sobolev space consisting of the set of square- integrable functions on S with square- intcgrablc 
weak derivatives of order smaller than or equal to k (see, for instance, Ref. [45| for more details). We define for each 
A > the following bounded, symmetric quadratic forms, 

Pa : H\S) x H l {S) -+ R, (/,«) ^ P x (f,u) :=< V/,Vu > +A < f,u >, 

Qa : H 2 (S) x P 2 (S) -> M, (/i, «) i— > Q A (ft, v) :=< (VV/i) TF , (VVw) TF > +A < V/i, Vw > +A < fc, w >, 

and the bounded linear functionals 

J w : H 1 (S) -> R, uH4(ti) :=<u,Vm>, 

K x : P 2 ^) -> R, v i-> A'x(w) :=< X, (Ww) TJ? >, 

With these notations, we reformulate the problems described in Eqs. (|A8IA9p as finding / <E H 1 (S) and h £ H 2 (S) 
such that 

P (/,u) = Ju,{u) for all u £ P 1 ^), 
Qo(M) = JfxW for all u e P 2 (5"). 

Notice that for A = 1, Pi is the standard scalar product on H 1 (S). Therefore, it follows by the Riesz representation 
lemma that the map L\ : H 1 (S) — > (H 1 (S))*, f H> P\(f, •) is invertible for A = 1. In particular, it is a Fredholm 
operator of index zero. On the other hand, by compact embedding, the difference operator L\ — Lq is compact, which 
implies that also the map Lq is Fredholm with index zero. Therefore, the equation Lof = J w has a solution if and 
only if Jcj(u) = for all u £ P 1 (5') such that V« = 0. Since this is always the case by the definition of J w , the 
existence of a solution / S H 1 (S) satisfying Eq. (|A8|) follows. By elliptic regularity theory (see, for instance Ref. (45| ) 
/ is C°° smooth and t\a '■= — Va/ is divergence-free. 

Similarly, Q\ is a scalar product on H 2 (S) for sufficiently large A > 0. In order to see this, we apply the integral 
identity (|A10|) below to (5a = V^/i, obtaining 

Q\(h, h) = \J (V A V B h)(V A V B h) + J (x - ^ (V A h)(V A h) + \J h 2 , 

s s s 

which shows positivity for A > fc/2 and also proves that in this case the norm induced by Q\ is equivalent to the 
standard norm on H 2 (S). Therefore, for A > fc/2 the maps M\ : H 2 (S) -> (H 2 (S))*, h M> Q\(h, ■) are invertible and 
Fredholm of index zero. Since M\ — Mo is compact, it follows that also Mq is Fredholm with index zero. Therefore, 
the equation M$h = Kx has a solution if and only if Kx(v) = for all v £ H 2 (S) such that (Wv) TF = 0, which is 
always the case by the definition of Kx- The existence of a solution h £ H 2 (S) satisfying Eq. (|A9j) follows, and by 
elliptic regularity it is C°° smooth. Then, Yab '■= Xab — (V Bh) TF satisfies the required property V a V b Yab = 0. 

This concludes the first part of the proof and shows existence and uniqueness of the decompositions ([A6IA7P . Notice 
that so far, no topological restrictions on S have been used. 

In the second part of the proof we show the existence of C°° functions g and k on S satisfying t]a = £A S V ' sg 
and Yab = £(A C ^ b)^ 'ck. In order to prove this, we apply the decomposition (|A6I) to the one-form a a ■= —eA B f]A- 
Therefore, a A = V A ff + Pa with \7 A /3 A = 0. Since = V A r/ A = V a e a b {V B g + Pb) = e AB V A f3 B , the one-form j3 B 
on S is both divergence- and curl-free. By the integral identity (|A10I) below, this implies (3 = since S has positive 
Gauss curvature, and it follows that t\a = £a B (Xb = £a S V sg- By a similar argument, we use the decomposition (|A7j) 
and write Yab = £(a C ? (V 'b)V 'ck + Z B ^ C ) where the C°° symmetric, trace-free tensor field Z satisfies V A V S Zab = 
and = S7 a \7 b Yab = V^V 5 'ea C Zbc- Therefore, the one-form r\A '■= V s Zab is both divergence- and curl-free 
which implies T] = 0. According to the second integral identity (|A11[) below, this implies Z = which concludes the 
second part of the proof. rj 

Lemma 1 Let (3 and Z denote a C°° one- form and symmetric, trace-free tensor field, respectively, on {S, g) . Then, 
the following integral identities hold: 



[{V A [3 B ){V A pB) + kp A f3 A ] 



\{V A p B - V B p A )(V A f3 B - V s /3a) + (V A f3 A ) 2 



(AlO) 



J [(V A Z BC )(V A Z B c) + 2kZ AB Z AB ] - 2 J g AB (V c Z CA )(V D Z DB ), (All) 



where k denotes the Gauss curvature of(S,g). In particular if k > it follows that there are no non-trivial divergence- 
and curl-free one-forms on (S, g) nor are there non-trivial symmetric, trace- and divergence-free tensor fields on (S, g). 
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Proof. For the first identity, we integrate 

(V A !3 B )(V A !3 B ) = (V A (3 B - V B (3 A )(V A (3 B ) + (V s /3 A )(V A p B ) 

over S and use Gauss' theorem twice to bring the second term on the right-hand side in the form of the square of the 
divergence of (3. For the second identity, we integrate 

(V A Z BC )(X7 A Z BC ) = (e AB V A Z BC )(e EF V E Z FC ) + {V B Z AC )(V A Z BC ) 

over S, use Gauss' theorem twice to rewrite the second term on the right-hand side as the square of the divergence of 
Z and use the fact that e A c Zc B = — Z A de d B ■ n 
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